/*  slfgcdl.cpp  by K.Tsuru  */
// function ID 2000   DRADIX, BRADIX
/***********************************************************************
SLong and SInteger classes
It provides the greatest common divisor by use of Knuth's binary method.
It returns a positive value.
See ReduceByPow2Rdx().
Ref: The Art of Computer Programing VOLUME 2 p. 338
************************************************************************/
#ifndef SN_H
#include "sn.h"
#endif
//sub-function : While "r" is an even number it divides by the power of two.
//The result is overwritten.
static void DivPow2(SLong& r){
    if(r[0] & 1u) return;
                                                    // 16                         8
    uint pow2 = (r.Type() == r.BIN_INT) ? (8u*sizeof(ulong)-BRADIX_BITS -1u) : 2u*DFIGURES;
    // = 16 for BRADIX, = 8 for DRADIX, 10000^2 = (2*5)^(2*DFIGURES)
    ulong div = 1uL << pow2;    //If the lowest two figures of "r" can be divided by this
                                //"r" can be too.
    while(div >= 2uL){
        while( !(r.Low2() % div) ) LDivPow2(r, r, pow2);
        div /= 2uL;     pow2--;
    }
#ifndef NDEBUG
assert(r[0] & 1);       // r is odd.
#endif
}
SLong gcdL(const SLong& x, const SLong& y){
    if( x.Type() != y.Type()) x.SetError(x.RADIX_ERR, "gcdL", 2000);
    SLong t(x.Type(), 0);
    if(x.IsOne() || y.IsOne()){     // x == 1 || y == 1
        t = 1.0;    return t;
    }
    SLong u(Labs(x)), v(Labs(y));

//a anterior process
//It divides by a common divisor 2^b*R^p to reduce the number of dividing by two
//in the next step.
    SLong divisor(x.Type(), 0);
    ReduceByPow2Rdx(u, v, &divisor);
#ifndef NDEBUG
assert( (u[0] & 1) || (v[0] & 1) );     // u or v must be odd.
#endif
    if(u[0] & 1u) t = -v;
    else          t = u;

    while( t.Sign(2000) ){
        DivPow2(t);     // while( !(t[0] % 2u) ) t = LDiv2(t);
        if(t.Sign() > 0) u = t;
        else             v = -t;
        t = u - v;
    }
//a posterior process
    if(!divisor.IsOne()) u *= divisor;
    return u;
}